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We investigate chemical bond formation and conductance in a molecular Cgo-junction under finite 
bias voltage using first-principles calculations based on density functional theory and nonequilib- 
rium Green's functions (DFT-NEGF). At the point of contact formation we identify a remark- 
ably strong coupling between the Ceo-motion and the molecular electronic structure. This is only 
seen for positive sample bias, although the conductance itself is not strongly polarity dependent. 
The nonequilibrium effect is traced back a sudden shift in the position of the voltage drop with 
a small Ceo-displacement. Combined with a vibrational heating mechanism we construct a model 
from our results that explain the polarity-dependent two-level conductance fluctuations observed 
in recent scanning tunneling microscopy (STM) experiments [N. Neel et al, Nano Lett. 11, 3593 
(2011)]. These findings highlight the significance of nonequilibrium effects in chemical bond forma- 
tion/breaking and in electron- vibration coupling in molecular electronics. 



PACS numbers: 73.63.-b, 68.37.Ef, 61.48.-c 

I. INTRODUCTION 

The influence of an external bias voltage and electronic 
currents on the formation and breaking of chemical bonds 
is a topic of increasing importance with the continued 
down-scaling of electronic components. This is especially 
accentuated in the limit of single-molecule devices.^ 

A substantial current may flow through a single bond 
and its effect on the stability and impact on transport is 
crucial. The phenomenon of random two- level conduc- 
tance fluctuations (TLF) is generally observed in a wide 
range of simple atomic and molecular contactsP^ It is 
often possible to relate these to changes in the bond- 
ing configuration driven by the current. Clearly, con- 
trolled and reversible switching between well-defined con- 
ductance states is a usef ul function.^ Over the years many 
examples of atomicPMD anc ] mo i ecu le-basedE3l2U switches 

have been demonstrated. However, the understanding 
of how the nonequilibrium electronic structure impact 
chemical bonding and conformational changes still pose 
many open questions. First-principles calculations and 
comparisons with well-characterized, time-resolved ex- 
periments can shed light on these issues. 

Nonequilibrium dyna mics of Ceo-systems has been un- 
der intense studjU^M^H. Here we focus on recently 
reported time-resolved measurements of single Cgo- 
contacts with a scanning tunneling microscope (STM)J2S1 
which showed that TLF occur in a narrow transition 
regime between tunneling and contact to Cgo- The ad- 
vantage of STM is the possibility to identify the orienta- 
tion of individual Cgo-molecules^ 2 ^ before and after con- 
trollable formation of the tip-molecule contact.^ More- 
over, the role of detailed electrode bonding geometry * 29 * 30 ! 
and contact point on the junction conductance has been 



clarified.^ 

More specifically, the experiment revealed the follow- 
ing interesting properties: (i) In the tunneling regime 
dl/dV spectroscopy shows that transport is dominated 
by the lowest unoccupied molecular orbital (LUMO) 
(seen as a resonance centered at a positive sample voltage 
of ~ 0.4 V), while (ii) in contact the I —V curve is close 
to linear in the voltage range [—0.4; +0.4] V, suggesting a 
relatively symmetric coupling of the LUMO resonance to 
the two electrodes. Intriguingly, (iii) the TLF was only 
observed at positive sample voltage around contact for- 
mation. These findings were discussed in Ref. 25 solely 
on the basis of dl/dV spectra in the tunnel regime. Es- 
sentially, only the spectral properties of the molecular 
adsorbate in equilibrium with the substrate were consid- 
ered. Here we present a different view on the experimen- 
tal findings based on our demonstration of a remarkably 
strong bias-dependent electronic coupling to the center- 
of-mass (CM) motion of the Cgo at the point when a 
bond is being formed between Cqq an d the apex atom 
of the STM tip. From first-principles calculations we 
obtain a detailed description of the Cgo-junction geome- 
try as well as the molecular LUMO resonance near the 
Fermi level. This allows us to construct a model for the 
TLF, which provides an explanation for the experimental 
findings. Our results demonstrate that the full nonequi- 
librium electronic structure needs to be accounted for to 
understand the observed TLF. 

Our paper is organized as follows. In Sec. [H] we de- 
scribe the first-principles method and our setup of the 
Ceo-contact system. In Sec. |III| we then describe the 
results obtained without fitting parameters for the con- 
tact formation between STM-tip and Cgo in equilibrium. 
Here we identify the formation of the chemical bond be- 
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FIG. 1: (color online) (a) Relaxed bond lengths, and (b) cor- 
responding total energy vs electrode separation L. The total 
energy is determined with respect to the initial configuration. 
The bond length di (tfe) between Ceo center-of-mass (CM) 
and surface (tip) and the junction geometry are defined in 
the inset in (a) along with the orientation of Ceo with respect 
to the underlying Cu surface (top right inset) and the tip apex 
atom (bottom right inset). The Ceo is oriented such that the 
tip apex atom approaches the 5:6 bond. 



tween the molecule and the tip apex atom. This is fol- 
lowed by our study of noncquilibrium effects and a dis- 
cussion of the identified polarity-dependent strong cou- 
pling between the Ceo CM motion and voltage drop 
(Sec. IV I. From these first-principles calculations we ex- 
tract in Sec. [V] parameters for a simple single-resonance 
model, most importantly the bias-dependent electron- 
vibration coupling to the CM motion. Together with a 
few additional parameters the model is used to calculate 
the TLF behavior, which can be compared to the exper- 
iment. Before concluding we discuss how the nonequi- 
librium forces modify the energy landscape for the CM- 
motion (Sec. VI). 



II. METHOD AND SETUP 

To study the contact formation and TLF we em- 
ploy the Siesta^ density functional theory (DFT) 
method, and its extension to finite bias using nonequi- 
librium Green's functions (DFT-NEGF) in the Tran- 
Siesta scheme !^ The generalized gradient approxima- 
tion (GGA-PBE) is applied for exchange and correlation 
(xc)P 



The Cgo junction geometry is modelled as shown in 
the inset to Fig. [T| The periodic supercell used in the 
DFT calculations contains one Cgo molecule supported 
on top of seven fixed Cu(lll) layers (27 Cu atoms per 
layer) with a pyramid-shaped Cu tip mounted on the bot- 
tom layer. To accurately describe the Cu surface and the 
chemical bonding with Cgo, an optimized diffuse basis 
set was applied for Cu surface layers and the tipp3 The 
counterpoise correction^ for the basis set superposition 
errors (BSSE) was applied to the total energy calcula- 
tions, which was checked against complementary calcu- 
lations with the Vasi^3 plane wave code as shown in 
Fig. 0b). 

The T-point approximation was employed for Brillouin 
zone integrations in the electronic structure calculation, 
while the transmission function was sampled over 3x3 
k-points in the 2D Brillouin zone parallel to the elec- 
trode surfaces. The residual atomic forces were lower 
than 0.02 eV/A for the atoms that were relaxed. The 
Ceo CM force constant was calculated from DFT total 
energies corresponding to configurations where the Ceo 
CM was rigidly displaced, up to 0.6 A from its equilib- 
rium position. 



III. CONTACT FORMATION 

We first focus on the bond-formation point at zero 
bias, and consider the approach of the STM tip towards 
a 5:6 Cgo-bond, i.e., a bond between a pentagon and a 
hexagon. We note that the fluctuations were observed for 
this orientation in the experiments and that no molec- 
ular rotations occur during contact formation in cither 
the experiments^ or j n our structure optimizations. 

We optimize the junction geometry by stepwise reduc- 
ing the size of the DFT-supercell in the direction per- 
pendicular to the surface, while relaxing the Cgo and tip 
atoms. Fig. [ija) shows the relaxed bond lengths d\ and 
d,2, between the Cgo center-of-mass (CM) and the surface 
and the tip atoms, respectively, function of elec- 

trode separation L. Around a characteristic separation 
Lq = 15.96 A, the distance di decreases rapidly while 
d% increases dramatically as the cell shrinks. This sig- 
nals the onset of a chemical bond formation between the 
STM tip and the Cgo molecule. This tip-Cgo attraction 
lowers the total energy of the system as witnessed by the 
binding energy curve in Fig. |T|(b) . 

The corresponding vibrational energy Ml associated 
with the Cgo CM motion as well as the zero-bias conduc- 
tance G =GqT(Ef) (conductance quantum Go = 2e 2 /h) 
of the junction are shown in Fig. [2] in the transition 
regime between tunneling and contact for the DFT equi- 
librium geometries. At L we find three eigenchannels 
contributing to the total transmission T{Ep) with the 
values {0.16,0.006,0.002}. The first channel dominates 
the transmission, because the three -fold degeneracy of 
the C 60 LUMO has been liftedPEH Thus, the C 60 sym- 
metry is broken in the contact configuration. One ob- 
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FIG. 2: (color online) (a) Vibrational energy of the Ceo 
center-of-mass (CM) motion between the electrodes, and (b) 
zero-bias conductance G of the junction vs electrode separa- 
tion L in the transition regime between tunneling and contact 
for the DFT equilibrium geometries. 



serves that the bond formation to the tip softens the Ceo- 
vibration [Fig. [2]ja)] and increases the conductance by 
roughly a factor of 2.5 [Fig. [2jjb)[. We note that the cal- 
culated conductance value of the order G = 0.2Go agrees 
very well with the experimental conductance in the tran- 
sition region between tunneling and contact where the 
TLF occurP^ Moreover, the calculated vibrational en- 
ergies agree with a recent theoretical study of the Ceo 
CM-motion on the Au(lll) surfaced 

According to our equilibrium DFT calculations we 
could not identify two well-defined stable configurations 
(for any fixed electrode separation) which could explain 
the existence of two different conductance states. Instead 
we observe a shallow energy landscape around the point 
of contact formation indicating that Cgo is rather free 
to move between the electrodes (e.g., the softening of 
the Ceo CM mode). We therefore speculate that a small 
barrier of the order of 10 meV, separating two distinct 
configurations, could be masked by limited numerical ac- 
curacy or by inherent approximations in the applied xc 
functional. In fact, recent theoretical studies of a some- 
what simpler system consisting of graphene on Ni(lll) 
have shown that various xc functionals can yield differ- 
ences in the potential energies describing the carbon- 
metal dist ance much beyond the energies relevant for 
our system.E!142l The disregard of current-induced forces 
acting on the atoms could also play an important role 
in the energy landscape, 43 a point we return to at the 
end of this paper. Finally, we note that the actual ex- 
periments involve a complex reconstructed surface struc- 
ture which we did not take into account. Because of 



these circumstances we shall therefore in our TLF-modcl 
(Sec. [v} postulate the existence of two configurations in 
the contact region separated by a small barrier (on the 
order of DFT-accuracy) , and instead focus our attention 
on the electron-CM vibration coupling and the resulting 
current-induced heating, which can explain the observed 
polarity-dependent TLF. 

As the electrode separation Lq is characteristic for the 
point of tip-Ceo bond formation, we take this configura- 
tion as the starting point for an exploration of how the 
nonequilibrium electronic structure and electron trans- 
port depend on Ceo motion. Fig. [3ja) shows the trans- 
mission spectra (with a prominent LUMO resonance) for 
several positions Adi of Ceo between the electrodes under 
three different applied sample voltages Vs- In each sit- 
uation the transmission fu nction is approximately given 
by a Breit-Wigner functiorPUil] 



T(E, V s ) « 2tt 



r r r 



T-L S 



I> + r. s 



[p T (E,V s ) + ps(E,V s )], (1) 



where 



Pt(s)(E,V s ) 



1 

2tt W- 



r T(g) 



•r s )/2] s 



(2) 



is the partial density of states of the LUMO resonance, 
positioned at eo(Vg), due to the coupling T T fg) to the 
tip (sample) electrode (neglecting energy dependence in 
Tt(s))- We take the equilibrium Fermi energy ep = 
as the energy reference and define the tip and surface 
chemical potentials as px — Vs/2 and ps = — Vs/2, re- 
spectively. The resonance parameters {£o,rT,Fs} are 
readily fitted to the DFT-NEGF calculations as a func 
tion of Cgo-position and voltage, as shown in Fig. |3 



IV. VOLTAGE DROP 

Remarkably, the nonequilibrium electronic structure 
reveals a strong variation of £o with Cgo-position for pos- 
itive sample voltages. This is a central finding of this 
work and below we shall show that it can explain the 
strong polarity dependence of the TLF seen in the ex- 
periments. In Fig. [3jc) we illustrate this by plotting the 
change in e relative to psi as a function of px for the 
various Cgo-displacements. For px < the eo mainly 
follows ps, while for px > a small increase in d\ and 
thus coupling to the tip, makes eo follow px rather than 
Ps, despite > Fx- 

The voltage dependence of eo, or equivalently the volt- 
age profile across the junction, can be understood roughly 
as a disposition of the system to maintain a constant elec- 
tron charge Q in the resonance.^ In order to illustrate 
this we consider a simple model calculation. Within the 
resonance model the LUMO charge is given by 



Q(Vs)= / PT (E,V s )dE 



p s (E,V s )dE. (3) 
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FIG. 3: (color online) (a) Transmission spectra T(E, Vs) 
for different sample bias voltage Vs as the Ceo-surface bond 
length is increased by Adi. Data points correspond to DFT- 
NEGF simulations and solid lines are fits to the Breit-Wigner 
formula [Eq. (TlJ>] for the transmission through a single molec- 
ular orbital, (b) Fitted LUMO level position eo an d coupling 
functions Ts and Ft due to the surface and tip coupling, 
respectively, and (c) voltage drop across the C6o-surface in- 
terface as a function of /it and Adi. The solid lines in (b-c) 
are guides to the eye. 
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FIG. 4: (color online) (a) Simple model calculation of the 
change in resonance position assuming charge-neutrality of 
the resonance at finite bias. We use so = 0.1 eV, Ts + Tt = 
0.2 eV (fixed), and vary Tt/Ts (distance to tip). The main 
voltage drop occurs between tip and Ceo (resonance) for Vs < 
and small tip coupling, while it shifts to the Ceo-substrate 
interface for Vs > and stronger tip coupling, (b) Voltage 
drop (change in DFT one-electron potential) calculated for 
increasing Ceo-surface bond length at a positive sample bias 
of 0.5 V. The voltage difference between consecutive contour 
lines is 18 mV (shown in a plane through the tip atom). 



If we assume constant LUMO charge independent of the 
applied bias, i.e., Q(Vs) — Q(0), we may determine the 
bias-dependent change in LUMO position, Seq(Vs) from 
Eq. ([3]). To mimic the change in Cgo-tip distance, d 2 , 
for fixed electrode distance, L. we vary Tt/Ts for fixed 
Tt + Tg. In Fig. [4^a) it is seen how this simplified model 
reproduces the cross-over in the full DFT calculation 
[Fig. [3|c)] for positive sample voltage when the contact 
is formed. Thus the main voltage drop changes from be- 
ing between tip and Ceo for Vs < to being between 
surface and C$o when Vs > and the distance to the tip 
is decreased (|fc — (is\ > |^ £ o — Mt|)- From Eq. ^ we 
can thus infer that in noncquilibrium there is a sensitive 
balance between coupling strengths (IV/s) and electrode 
chemical potentials (ht/s) that can displace the voltage 
drop from one interface to the other with a small relative 
change in coupling strengths. 

The voltage drop effect can also be seen directly in 
the actual voltage drop landscape (change in the one- 
electron potential with respect to equilibrium) shown in 
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Fig. |4jb). The voltage drop is observed to shift from the 
Cgo-tip interface to the Cgo-substrate interface with a 
small Ceo-displacement, an effect not present for Vs < 
(not shown). 



V. HEATING AND FLUCTUATIONS 

We next explain how the strong variation of £q with 
Cgo-position for Vs > can be related to the strong 
polarity dependence of the TLF. We start by assuming 
that the main current-dependence comes from the exci- 
tation of Cgo CM-motion, described by a harmonic po- 
tential with Ml ss 4 meV [cf. Fig. |2ja) at Lq\. Guided 
by the fact that the switching rates observed in the ex- 
periments (ms time scale) are very slow compared to CM 
oscillations, we propose that the switching involves a slow 
"bottle-neck" process, possibly involving tunneling along 
the reaction coordinate (RC), and that this process takes 
place when the excursion of the Cgo (Ae^) is beyond some 
critical distance from the equilibrium position. Inspired 
by the study of tunneling of a Cgo molecule in the low- 
conductance regime^ we express the switching rate as 



R(Ys) 



-dl/iA^iVs) 2 ) = 



-A/(Bcm(Vs)) 



(4) 



i.e., as a product of the probability of Cgo being at an 
excursion Adi = d c away from equilibrium and of a rate, 
r s , describing the slow process along the RC. The crit- 
ical distance d c , or equivalently the energy barrier A, 
controls how far the Cgo needs to move in order to facili- 
tate switching. The mean displacement (Adf), or equiv- 
alently the mean oscillator energy (Ecm(Vs)) , are quan- 
tities which we can calculate within our TLF-modcl. 

The excitation of the Cgo-CM motion by the current is 
determined from the electronic coupling to this motion. 
Using Fig. |3jb) we extract the electron-vibration cou- 
pling M from the shift in resonance position with Ceo- 
displacement d\ via^ 



M{V S 



V2 



(3) 



We evaluate the slope, d^SoiVs), around Arfi = 0.3 A, 
which is in the middle of the transition region [Fig.[3jh)], 
and note that the slope does not change significantly as 
we inc rease A d±. The characteristic oscillator length is 
la = y/h/mVl 0.04 A (Cgo mass m), which is compara- 
ble to the size of the transition region in Fig. [2j The ex- 
tracted electron- vibration coupling, M, is shown in Fig. [5] 
as a function of sample voltage. A remarkably strong en- 
hancement is evident for Vs > 0. 

The excitation of the CM-motion, as seen in its 
mean energy (Ecm(Vs)) , can be obtained from the bias- 
dependent rates of phonon emission, 7cm(V's'), and of 
electron-hole pair generation, 7 h(Ks')- These rates can 
be determined within first order perturbation theory 



(Fermi's Golden rule). Since Ml is much smaller than 
all other electronic parameters, we may write 



7em(V S ) 



4tt 



\M(V S )\ 2 0( 



e\Vs\ 



i 



(G) 



r\V s \/2 

I P s(E,Vs)p T {E,Vs)dE, 

-|Vs|/2 

|2 , 



7ch(^s) « 4ir\M(V s )[ z n (7) 
x [p(ht, V s )pt(pt, V s ) + p(ps, V s )ps(ps, Vs)] , 

where p = p? + ps- From these rates we can write a rate 
equation for the mean phonon occupation (n(Vs)), 

(n(V s )) =7em(Vs)-{7«o»(Vs)+7ph}{<n(V s ))-n B } > (8) 

where 7 p h represents the vibrational relaxation due to 
anharmonic coupling to phonons in tip / substrate and ns 
is the Bose-Einstein (equilibrium) phonon occupation of 
the considered mode. The steady-state solution is simply 



(n(V s )) = n B + 



7em(Vs 



7ch(^s) +7ph' 



(9) 



Following Refs. 47 48 one can estimate a phonon damping 
to the substrate of Cgo-CM motion via the formula 



7 P h = — |^ 4 «0.1f7 
m Cu 2oj% 



(10) 



where mcu is the mass of a substrate atom, and w e = 30 
meV a frequency characterizing the elastic response. This 
damping is likely to be exaggerated compared to the ex- 
perimental situation since the Cgo is adsorbed on a re- 
constructed surface with low-coordinated surface atoms 
and lower density of long wavelength phonons. This is 
a critical point for the explanation of the experimental 
result. We find that the best agreement is obtained for 
7ph ~ 0.00157. In Fig. [6] we show how 7 cm , 7ch, and (n) 
varies with the sample voltage Vs along with the effec- 
tive temperature defined through a Bose-Einstein distri- 
bution (n(V s )) = l/{ e hn / k ^ T ^ v s) - 1). In all cases we 
see an enhancement for Vs > 0. If we use 7 p h sa 0.1 Q 
the mean occupation and effective temperature become 
a factor 100 smaller, but exhibit the same behavior as 
in Figs. [6jc)-(d). Finally we can calculate the oscillator 
energy as {E C m(V s )) = Ml((n(V s )} + 1/2) and thus the 
current-dependent rate from Eq. Q. 

Fig. |7ja) shows how the calculated switching yield 
Y(Vs) (blue squares), defined as the switching rate per 
tunneling electron, can reproduce the experimental data 
(red circles) if we use a barrier height of A = 2MI, a 
"tunnel-rate" r s — 80 ms^ 1 , 7 p h ~ 0.0010, and a back- 
ground temperature of T = 14 K as fitting parameters. 
The slightly elevated temperature, compared to the ex- 
periment performed at T — 7 K, helps to smoothen the 
onset of the switching rate at small voltages. This can be 
justified by vibrational heating of other modes and their 
anharmonic coupling to the CM motion of the Cgo- The 
relatively slow r s is consistent with a tunneling process, 
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FIG. 5: (color online) Calculated electron-vibration coupling 
\M\ for the Ceo CM motion (black squares) as a function of 
sample voltage Vs, based on fits to the LUMO resonance in 
the transition region at the electrode separation Lq. 
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FIG. 6: (color online) Rates of (a) electron-hole pair gener- 
ation , (b) phonon emission, and resulting (c) mean phonon 
occupation as a function of sample bias for 7 P h ~ 0.00117. 
(d) Effective temperature T e g[Vs) corresponding to the mean 
phonon occupation in (c). 



FIG. 7: (color online) (a) Calculated switching yield Y(V$) 
(blue squares) based on the excitation of the Ceo CM-motion 
[Ml = 4 meV) using Eq. p} with r s = 80 ms~\ T = 14 K, 
and A = 2Ml. (b) Switching yield in the case of constant 
electron-vibration coupling, M[Vs) = M[Vs = 0). (c) Same 
as (a) but with 7 P h ~ 0.117 <C 7 c h. For comparison also the 
experimental data from Ref. 1251 are shown (red circles). 



and the small A is consistent with the fact that we could 
not determine the barrier with our DFT calculations. 

In Fig. [7^b) we show the calculated fluctuation rate 
in the case of constant zero-bias electron-vibration cou- 
pling where only the spectral energy-dependence of the 
molecule are considered, cf. the explanation presented 
in Ref. H5] for the polarity dependence. However, it 
is clear that we are only able to reproduce the experi- 
mental results if we take the behaviour of the electron- 
vibration matrix element with bias into account. These 
findings suggests that (i) the strong polarity dependence 
of the switching is rooted in the noncquilibrium electron- 
vibration coupling in the transition region where the 
bond-formation between tip and Cgo takes place, and 
(ii) that the reason for the observed saturation of the 
switching rate per electron is due to the steadily increas- 
ing electron-hole pair damping with bias, Eq. (j8j, so 
this becomes comparable with 7 p h- This is an impor- 



tant point as illustrated in Fig. 7(c) where we show how 
the switching yield using the estimate 7 p h = 0.1f2 grows 
for Vs > 0.2 V (contrary to the experiment). 

We note that the calculated current is roughly linear in 
voltage as in the experiments, and thus do not contribute 
significantly to the polarity dependence of the switching 
compared to the pronounced effect seen in Fig. [5] for the 
electron-vibration coupling M. We further note that one 
theoretical stud}®' has previously reported a nonlinear, 
polarity- independent M for a smaller symmetric molec- 
ular junction and only at significantly higher voltages 
V > 0.4 V. 
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VI. EFFECT OF CURRENT-INDUCED FORCE 

In this section we estimate the change in the poten- 
tial energy landscape in Fig.JIJb) when a nonequilibrium 
force is exerted on the Cgo~tip bond during contact for- 
mation. To calculate this additional bond force we con- 
sider the interaction between the Cgo LUMO resonance 
and a wide band centered on the tip (see inset in Fig. [8]) . 
For this system we define the following two-site Hamilto- 
nian 



H 



e (d 2 ,V s ) t(d 2 ,Vs) 
t(d 2 ,V s ) MVs) 



(11) 



where we explicitly stress the dependence on bias and 
bond length d 2 . The interaction strength t(d 2 , Vs) is cal- 
culated using r^(c?2, Vs) = 27T7y|i(<i2, Ks)| 2 , where jt is 
the wide band density of states on the tip, i.e., a con- 
stant which can be fitted to reproduce the transmission 
spectra in Fig. ^a). The bond force is calculated using 
the general expression^! 

F b (d 2 ,V s ) = -2Tr[(d rf2 H)D] (12) 

= -2 [Cii9 d2 e (d 2 , V s ) + 2D l2 d d2 t{d 2 ,V s )] , 

where a factor of 2 is included to account for spin. The 
elements, Dn and D\ 2 {= D 2 i) of the density matrix D 
are determined from the spectral properties of the con- 
sidered states,^ which can be calculated from the fits in 
Fig. jib). 

Since we only consider motion along a single coor- 
dinate the current-induced force is energy conserving, 
SF b (d 2 ,V s ) = F b (d 2 ,V s ) - F b (d 2 ,0), enabling us to cal- 
culate the change in bond energy at a given bias, 



AE b (d 2 ,V s ) = - SF b (z,V s )dz. 



(13) 



The integration limits are defined such that initially a 
contact is established, d 2 ^ ~ 5.6 A cf. the equilibrium 
curve in Fig.[8j and the n we integrate along a path where 
the Ceo-tip contact is gradually separated. Addition of 
the energy term in Eq. ( |13| on top of the equilibrium 
total energy in Fig. [ljb) yields the bias-corrected curves 
shown in Fig. [8] Astoundingly, we see that only at posi- 
tive Vs a tiny barrier of the order a few meV may appear 
between two stable configurations corresponding to con- 
tact and tunneling cases, respectively. The origin of the 
significant lowering of the tunneling part of the binding 
energy curve for positive Vs is related to the asymmetry 
in resonance position, which yields a large contribution 
from dd 2 eo in Eq. ( 12 1 only at positive Vs. Finally, we 



note that the order of magnitude of the nonequilibrium 
barrier is in accordance with our assumption in the fluc- 
tuation calculation in Fig. [7] 



VII. CONCLUSIONS 

In summary, we have presented the results of first- 
principles calculations which combined with a heating 




5.5 6.0 6.5 7.0 7.5 8.0 

cfe (A) 



FIG. 8: (color online) Equilibrium DFT total energy from 
Fig. [ljb) plus an additional contribution due to a (conserva- 
tive) current- induced force. Note that energy is here plotted 
against Ceo-tip bond length di instead of electrode separation. 
Inset: Simple model describing a bond between the molecular 
resonance at eo and a wide band centered on the tip, which 
follows the chemical potential /it- 



model and assuming a small energy barrier can explain 
the experimentally observed bias-dependent TLF ob- 
served in a Ceo STM junction. Our main point is that 
the electron- vibration coupling can depend very strongly 
on the bias polarity. In this system we can trace this 
back to sensitivity of the nonequilibrium electronic struc- 
ture/voltage drop with respect to the Cgo-motion just 
when the contact is being formed. The bias dependence 
of the electron- vibration coupling has so far not been con- 
sidered in most calculations of inelastic electron trans- 
port and current-induced excitations. It remains to be 
answered to what extend this is important in general. In 
order to model the experimental switching we had to as- 
sume a small energy barrier for the Ceo-motion at the 
contact formation point. Although it is likely that the 
small barrier is masked by inaccuracy inherent in the 
DFT, the finite unit-cell employed, or numerical error, 
we showed that the nonequilibrium can induce significant 
changes in the potential energy surface. Our estimate of 
the current-induced force exerted on the Cgo-tip bond 
did indeed indicate an energy barrier for positive sample 
voltage. 

In the presence of a significant current a number of 
different excitation mechanisms can become active. Re- 
cently, it has been discussed how current-induced forces 
can lead to "run-away" instabilities such as bond-rupture 
for highly conducting systems G ~ Go, and volta ges in 
the range involved in the present experiment! 52 ! 53 ! . TLF 
experiments seems to be a promising way to probe these. 
The runaway effect requires the action of several vibra- 
tion modes and we have limited our discussion here to 
a single mode. Our results demonstrate how the full 
nonequilibrium electronic structure can be of crucial im- 
portance for the formation/breaking of chemical bonds 
and electron-vibration coupling in the presence of cur- 
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rent. 
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